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Abstract 

The crucial role that Ions Energy Distribution Function (IEDF) at 
the electrodes plays in plasma processing of semiconductor materials de- 
mands that this quantity be predicted with high accuracy and with low 
noise levels in any plasma simulator. In this work, an efficient Particle- 
in-cell/Monte-Carlo (PIC/MC) simulator is developed to model IEDF at 
the electrodes of electronegative plasma discharges. The simulator uses an 
effective method to increase the number of MC particles in regions of low 
particle density by splitting the particles and by adjusting their statistical 
weight. This statistical enhancement technique, which does not require 
interprocessor communication, is particularly suitable for parallel process- 
ing. The simulator is used to model an electronegative rf discharge at a 
pressure of 25 mTorr. The IEDF obtained from this simulator has good 
statistics with low noise levels, whereas the IEDF calculated by standard 
PIC/MC simulator is jammed with stochastic noise. 



1 



I. Introduction 



Radio frequency glow discharges are widely used in microelectronics for de- 
position and etching thin solid films and for Ultra Large Scale Integrated (ULSI) 
circuit fabrication |Q . The understanding of these discharges is of considerable 
importance for a better control and optimization of semiconductor processing. 
Unfortunately, these discharges are complex systems that are very difficult to 
analyze [2-4]. In the last few years, Particle-In-Cell/Monte-Carlo (PIC/MC) 
simulation has been used successfully to model rf discharges. PIC/MC tech- 
nique combines the microscopic transport of electrons and ions with the self 
consistent electric field. The technique is very attractive because it can handle 
the nonlocal effects which are dominant at low pressures without making any 
ad hoc assumptions [5-6]. Furthermore, the PIC/MC technique is capable of 
predicting ions energy distribution function at electrode surfaces. 

Unfortunately, PIC/MC simulation of RIE discharges suffers from signifi- 
cant noise which is inherently associated with regions of low particle density. 
In particular, because the particle density in the sheath is many orders of mag- 
nitude smaller than in the bulk, the MC noise inevitably associated with the 
sheath is considerably larger than the noise associated with the bulk. Several 
techniques have been proposed to improve the sampling within the sheath [7-9] . 
However, even with the improved sampling techniques, the present PIC/MC 
simulators suffer from intensive computational requirements and from excessive 
noise associated with the sheath. Serious problems arise from the noise within 
the sheath and they may threaten the reliability of the simulation data. 

II. Details of simulation 

Since the motion of electrons and ions in rf discharges under the influence of 
self-consistent electric field bears some resemblance to the motion of electrons 
and holes in semiconductor devices, the advanced techniques that have been 
developed for device simulation [10-12] may also be useful for plasma simulation. 
For example, for hot electron injection in MOSFETs (r|, only one out of 10 10 
channel electrons makes a contribution to the gate current. A straight forward 
MC simulation is extremely expensive, requiring the simulation of billions of 
cold MC particles. To reduce the simulation cost of these problems, Phillips 
and Price JhJ divided the domain of the electron states S into a rare part R 
and a complementary common part C. When a particle enters R from C, its 
state is stored. The rare portion of the particle history is then repeated for 
a fixed number of N times starting from the same initial condition. After N 
repetitions, the MC simulation continues in the C domain using only one of the 
iV rare trajectories. Obviously, in calculating the average values, a weight of 
1 /N is given to the results in R. 

Motivated by their work, we use a similar technique to increase the number 
of MC particles in the sheath by splitting the particles in regions of low particle 
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density. There are, however, some fundamental differences between this work 
and that of [10-12]. In particular, in the work of Phillips and Price, the energy 
components of the particles are weighed differently, whereas in the present work, 
as is shown in part III, the spatial components of the particles are weighed 
differently. 

There are also some important differences between this work and the previ- 
ous plasma simulation works that increase the number of MC particles within 
the sheath [14-16]. For example, in the work of Goto et. at the dis- 

charge is divided into many slabs and every slab consists of M s sub-slabs. Each 
sub-slab in turn is divided into Mf fine-sub-slabs in which the generation and 
loss of particles are counted. The number of particles in every fine-sub-slab is 
set to a standard number N s . When the simulated particle number of a slab 
exceeds 1.4M s MfN s , the number of particles of all fine-sub-slabs in the slab is 
counted. If the particle number of a fine-sub-slab is greater than N s , then it 
is reduced to N s . Thus in their work, as well as in Refs. [14-15], an attempt 
is made to set the number of particles in every fine-sub-slab or in every cell 
equal to a constant number N s . In contrast, in this work, the discharge is not 
divided into slabs, sub-slabs, and fine-sub-slabs, and no attempt is made to set 
the number of particles in every cell or in every fine-sub-slab equal to a standard 
number. In the present work, a suitable weight Wi is associated with each cell i 
which is automatically updated with every time step, and unlike Refs. [14-16], 
the particles are split when the particle's weight is significantly larger than the 
weight associated with the cell. In this respect the present work is similar to 
the splitting technique used in the device simulation [10-12]. 

III. Simulation of electrongative discharge with 
SPLIT-PIC/MC technique 

A question of fundamental interest is whether the SPLIT technique, which 
was applied to electropositive discharges [0, can also be used to improve the 
statistics of electronegative simulators. To address this question, we developed 
a simulator to model electronegative gases [18-20] using SPLIT technique. The 
present simulator uses an effective method to increase of the number of MC 
particles in regions of low particle density by splitting the particles and by 
adjusting their statistical weight. The computation time for splitting the MC 
particles is very small and in general is less than 5% of the total calculation 
time. Using this splitting technique, the number of MC particles within the 
sheath becomes of the same order as in the bulk. In contrast, in the standard 
PIC/MC simulation, even though the sheath is the most important region in 
the discharge and has a profound influence on the properties of the bounded 
plasma, the number of MC particles in the sheath is many orders of magnitude 
smaller than in the bulk. 

Very briefly, the present SPLIT- PIC/MC simulator consists of the following 
steps: 
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(1) Each cell i contains a suitable weight Wi which is automatically updated 
with time. 

(2) The instantaneous locations of MC particles representing ions and electrons 
are interpolated to the grid points to obtain the charge density. 

(3) Poisson's equation on a spatially discretized mesh is solved to obtain the 
electric field. 

(4) The electric field is interpolated from the grid points to the location of MC 
particles. 

(5) The equations of motion under the local and instantaneous electric field are 
integrated. 

(6) When ion j enters cell i, its weight rrij is compared with Wi. If the particle's 
weight is much heavier than Wi, then it splits into two sub-particlcs, each with 
weight 

(7) Random numbers (Monte Carlos technique) and collision cross sections are 
used to determine the probability that each particle suffers a collision. 

It is worth noting that the present SPLIT-PIC/MC simulator is particu- 
larly suitable for parallelization and vectorization. The SPLIT method does 
not use the distribution of particle weights in a cell to adjust the number of MC 
particles to a desired population. Since the SPLIT technique does not require 
fetching particle data handled by other processors, there is no need for inter- 
processor communication, and complete parallelization with respect to particles 
is achieved. The numerical results presented here were all performed on an 
NEC parallel machine Cenju-3 (16 processor) which has a VR4400SC (75 MHz) 
processor and 64MB local memory in each processor [^lj . 

To examine the effectiveness of SPLIT technique, we have modeled an elec- 
tronegative rf plasma at a pressure of 25 mTorr both with the SPLIT-PIC/MC 
simulator and with the standard PIC/MC simulator p2| . The simulator is 
based on a Chlorine- like gas. The processes considered are ionization leading to 
positive ion formation, dissociative electron attachment leading to negative ion 
formation, and positive ion/negative ion recombination leading to positive and 
negative ions removal. These processes are represented as follows: 

Cl 2 + e C7+ + 2e (1) 
Cl 2 + e -> Cr + CI (2) 

cii + cr -► ci 2 + ci (3) 

where C7^~ and Cl~ are the main positive and negative ions respectively. In the 
simulation, the left electrode is driven with a voltage V r f(t) = Vrfsinut, where 
V r f = 300 Volts, the applied frequency is -^^r = 10 MHz, and the discharge 
length is 10 cm with perfectly absorbing electrodes. The time step St — 0.1ns 
which is small enough to resolve the electron plasma frequency. 

Since in this work we are primarily interested in testing the effectiveness 
and usefulness of SPLIT-PIC/MC simulator, we have intentionally kept the 
model of the plasma simple. The electron-neutral processes consist of elastic 
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scattering and ionization. The total electron-neutral scattering cross section 
o 'total (v ) is given by <7 to tai (v ) = ^ t ° tal , where w is the electron velocity and the 
rate constant K to tai is 5 x 10 -8 cm 3 / s ||. Ionizing collisions occur if the elec- 
tron energy is larger than 16 eV. An ionizing collision [Eq. (1)] is modeled by 
loading a new electron and ion at the position of the ionizing electron. The ki- 
netic energy after ionizing collision is partitioned between the two electrons with 
equal probability. The cross section for ion-neutral charge-exchange collision is 
c C e(w) = 1 x 10~ 15 cm 2 ||. After charge exchange collision, the ion is scattered 
isotropically, with its energy being set at the background temperature. Since the 
rate constant for attachment K at is constant |Q, the cross section for attach- 
ment [Eq. (2)] is modeled as a at (v) = ^f£ where K at = 1.8 x 1CT 10 cm 3 /s @. 
The attachment reaction is modeled by removing the electron from the system 
and by loading a negative ion at the position of the attaching electron. The ki- 
netic energy of the generated negative ion is set at the background temperature. 
Finally because the rate constant for recombination K at is also constant, the 
cross section for recombination [Eq. (3)] is modeled as <r rec {v) = ^ rec , where 
Krec = 5 x 10~ 8 err? I s Jig] , and Vi on is the velocity of the recombining ion. The 
recombination reaction is modeled by removing positive and negative ions from 
the system. Simulation results have shown that the splitting technique leads to 
reduction of noise within the sheath no matter what the detailed forms of the 
cross sections are. 

IV. Results and discussion 

The simulation is started by loading a few hundred quiet particles uniformly 
between the electrodes. Initially the velocity of all particles, negative ions, 
positive ions and electrons, is set to zero and the density of positive ions is set 
equal to the density of negative ions and electrons. Since the mass of electrons 
is much smaller than the mass of ions, electrons rapidly leave the system and 
sheaths develop at the boundaries. The density of electrons and ions continues 
to rise by ionization until steady state is reached where the ionization rate 
becomes equal to the total loss rate. At steady state the number of electrons 
and ions approaches a constant value. With SPLIT-PIC/MC technique, the 
simulator approaches steady state after approximately 900 cycles, whereas with 
the standard PIC/MC simulator, more than 2000 cycles are required for the 
simulator to arrive at steady state. 

To clearly demonstrate the effectiveness of the SPLIT technique, we ran 
SPLIT-PIC/MC simulator for approximately 1000 cycles, while we ran the stan- 
dard PIC/MC simulator for approximately 3000 cycles. Since the process of 
splitting particles takes less than 5% of the calculation time, the total computa- 
tion time for the SPLIT simulator is approximately three times smaller than for 
the standard simulator. As shown below, even with 66% smaller computation 
time, the simulation results obtained from SPLIT-PIC/MC simulator are more 
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accurate and have lower noise levels than the results obtained from the standard 
PIC/MC simulator. 

Figures 1 (a) and 1 (b) show the number of MC particles for positive ions 
obtained from the SPLIT-PIC/MC simulator and from the standard PIC/MC 
simulator (in relative units). The important point to notice is that for SPLIT- 
PIC/MC simulator, the number of MC particles in the sheath is of the same 
order of magnitude as in the bulk. In contrast, for the standard PIC/MC sim- 
ulator, the number of MC particles in the sheath is many orders of magnitudes 
smaller than in the bulk. 

Figures 2 (a) and 2 (b) show the instantaneous (not averaged over several 
period) positive ion density in log scale obtained from the SPLIT-PIC/MC 
simulator and from the standard PIC/MC simulator. Note that the massive 
ions do not respond to the rapid variations of the electric field and hence the 
ion density profile is nearly time independent. Figs. 2 (a) and 2 (b) clearly 
indicate that SPLIT technique leads to significant reduction of MC noise within 
the sheath. Figure 2 (c) shows the positive ion density profile within the sheath 
in linear scale simulated by the standard PIC/MC technique (dashed) and by the 
SPLIT-PIC /MC technique (line) . Comparing the two cases, it can be seen that 
the PIC/MC profile is zero within several cells, whereas the SPLIT-PIC/MC 
profile is never zero within any cells and has a much smaller noise level. 

Figure 3 shows the instantaneous profile of negative ion density within the 
discharge simulated by SPLIT-PIC/MC technique. The negative ions are con- 
fined to the bulk of the discharge and are completely excluded from the sheath 
by the large sheath electric field. 

Figures 4 (a) and 4 (b) show the instantaneous positive ion phase-space sim- 
ulated by the SPLIT-PIC/MC technique and by the standard PIC/MC tech- 
nique. These figures clearly show that the ions are accelerated by the sheath 
electric field. The ions that suffer several collisions arrive at the electrode with 
very small energies. In contrast, the ions that do not suffer any collisions arrive 
at the electrode by the average sheath potential. The important point about 
these figures is that the SPLIT-PIC/MC simulator dramatically increases the 
number of ions within the sheath, hence leading to significant reduction of MC 
noise in regions of low particle density. 

As another example of the effectiveness of SPLIT technique, the variation 
of the electric field within the discharge at phase 4^ is simulated both by the 
SPLIT-PIC/MC [Fig. 5 (a)] technique and by the standard PIC/MC [Fig. 5 
(b)] technique. Due to disparity between electrons and ions masses, electrons 
rapidly leave the system and large electric fields develop near the electrodes. 
The electric field in the bulk of the discharge, however, is much smaller than in 
the sheath, indicating that the net charge density within the bulk is very small. 
The large electric field at the electrode boundaries effectively pushes the negative 
ions out of the sheath (see Fig. 3). Figure 5 (a) clearly demonstrates that 
SPLIT-PIC/MC technique, by increasing the number of MC particles within 
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the sheath, dramatically reduces the noise associated with the sheath electric 
field [it should be emphasized that these figures are not averaged over many 
periods, rather they represent the snapshot or the instantaneous electric field. 
Furthermore, the computation time for Fig. 5 (b) is three times larger than for 
Fig. 5 (a)]. 

To demonstrate that the present SPLIT-PIC /MC simulator can predict the 
physical properties of the electronegative discharges, in Figs. 6 (a) and 6 (b), 
we present positive and negative ion current densities averaged over 50 periods 
and at four times during the rf cycle. (In contrast to Figs. 2-5, Figs. 6(a) and 6 
(b) are averaged over many periods and are not instantaneous). Note that the 
negative ion current density is almost one order of magnitude smaller than the 
positive ion current density. Furthermore, since negative ions are excluded from 
the sheath, the negative ion flux within the sheath is almost zero. In contrast, 
positive ions are accelerated by the sheath electric field and positive ion flux 
increases sharply within the sheath. 

Figures 2-5 have clearly shown the effectiveness of SPLIT-PIC/MC simula- 
tor for improving the statistics of ion density and the discharge electric field. 
However, from the point of view of semiconductor processing, the most impor- 
tant characteristic of RIE discharges is the Ions Energy Distribution Function 
(IEDF) at the surface of the electrodes [23, 24]. Thus a question of fundamental 
interest is whether the SPLIT technique can improve the statistics of IEDF at 
the electrodes (the question of the effectiveness of SPLIT technique on IEDFs 
at the electrodes was not considered in |l7|]). To answer this question, we sim- 
ulated the IEDF both by the standard technique and by the SPLIT technique. 
The IEDF obtained from the SPLIT-PIC/MC simulator and from the standard 
PIC/MC simulator are shown in Figs 7 (a) and 7 (b). Figure 7 (a) is obtained 
after running SPLIT-PIC/MC simulator for 950 periods, whereas Fig. 7 (b) is 
calculated after running PIC/MC simulator for 3000 periods; however, in both 
Figs. 7 (a) and 7 (b), the IEDFs are time averaged over 50 periods. The IEDF at 
the electrode consists of a few spikes that are superimposed on a collisional back- 
ground. Although most ions have many collisions in the sheath and hence create 
the collisional background profile, a few ions have only several charge-exchange 
collisions, creating the spikes in the profile. Comparing Figs. 7 (a) and 7 (b), it 
can be seen that with the standard PIC/MC simulator, the statistics of IEDF 
is very poor, jammed with stochastic noise; whereas with the SPLIT-PIC/MC 
simulator, the statistics of IEDF is significantly enhanced and has a much lower 
noise level. To further demonstrate the effectiveness of the SPLIT technique, 
in Figs. 7 (c) and 7 (d), we present the IEDF at the electrode simulated by 
the standard PIC/MC technique and time averaged over 150 and 600 periods 
respectively. Note that Fig. 7 (a) which is obtained from SPLIT-PIC/MC sim- 
ulator and is time averaged over only 50 period is smoother than Fig. 7 (d) 
which is obtained from the standard PIC/MC simulator and is time averaged 
over 600 periods. These figures clearly demonstrate that SPLIT technique leads 
to dramatic reduction of stochastic noise in IEDF at the electrodes. 
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V. Conclusion 



A stable algorithm which is highly suitable for parallelization and vectoriza- 
tion is developed to simulate electronegative rf discharges. The algorithm splits 
Monte-Carlo particles in regions of low particle density and hence improves the 
sampling and reduces the noise inherently associated with the sheath. The pro- 
cess of splitting the particles takes only a very small fraction (less than 5%) of 
the total computation time. The most attractive feature of SPLIT technique is 
that it can dramatically reduce the noise associated with the IEDF at the elec- 
trode. In particular, the IEDF at the electrode obtained from SPLIT-PIC/MC 
simulator and averaged over 50 periods has a better statistics and lower noise 
level than the distribution function obtained from the standard PIC /MC simula- 
tor and averaged over 600 periods. The SPLIT-PIC/MC simulator shows great 
promise for modeling any kind of bounded plasmas where sheaths inevitably 
develop, and for calculating IEDF at the electrodes of more complex discharges 
(such as ECR or ICP) in two or three dimensions. 
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Figure Captions 

Important Note: All figures should be considered from left to right and top 
to bottom. 

Figures 1 (a), 1 (b), are on the same page. 

Figures 2 (a), 2 (b), and 2 (c) are on the same page. 

Figure 3 is on one page. 

Figures 4 (a) and 4 (b) are on the same page. 
Figures 5 (a) and 5 (b) are on the same page. 
Figures 6 (a) and 6 (b) are on the same page. 
Figures 7 (a) and 7 (b) are on the same page. 
Figures 7 (c) and 7 (d) are on the same page. 

Fig. 1 (a) - Number of MC particles for ions obtained from SPLIT-PIC/MC 
simulator. Simulation conditions: V r f — 300 Volts, = 10, pressure 25 mTorr. 
These conditions apply to all figures. In Figs. [1-6], x — (left electrode) is the 
powered electrode and x — 10 (right electrode) is the grounded electrode. 
Fig. 1 (b) - Number of MC particles for ions obtained from the standard 
PIC/MC simulator. 

Fig. 2 (a) - Instantaneous ion density obtained from SPLIT-PIC/MC simulator 
in log scale. 

Fig. 2 (b) - Instantaneous ion density obtained from the standard PIC/MC 
simulator in log scale. 

Fig. 2 (c) - Instantaneous ion density in the sheath from PIC/MC simulator 
(dashed) and from SPLIT-PIC/MC simulator (line). 

Fig. 3 - Instantaneous negative ion density within the discharge obtained from 
SPLIT-PIC/MC technique in log scale. 

Fig. 4 (a) - Ion phase space diagram obtained from SPLIT-PIC/MC simulator. 
Fig. 4 (b) - Ion phase space diagram obtained from the standard PIC/MC 
simulator. 

Fig. 5 (a) - Instantaneous electric field at phase ^ obtained from SPLIT- 
PIC/MC simulator. 

Fig. 5 (b) - Instantaneous electric field at phase ^ obtained from PIC/MC 
simulator. 

Fig. 6 (a) - Positive ion current density averaged over 50 periods at 4 times 
during the rf cycle. For Figs. 6 (a) and 6 (b) ^ = (solid line), ^ = 0.25 
(dotted line), ^ = 0.5 (dashed-dashed line), ^ = 0.75 (dashed line). 
Fig. 6 (b) - Same as Fig. 6 (a) but for negative ion current density. 
Fig. 7 (a) - Ions energy distribuion profile at the electrode obtained from 
SPLIT-PIC/MC simulator after running the simulator for 1000 periods and 
averaging over 50 periods. In Figs. [7 (a), 7 (b), 7 (c), and 7 (d)], the horizontal 
axis represents the energy of ion striking the electrode, and the vertical axis 
represents the distribution of ions energy. 

Fig. 7 (b) - Ions energy distribuion profile at the electrode obtained from 
PIC/MC simulator after running the simulator for 3000 periods and averaging 
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over 50 periods. 

Fig. 7 (c) - Ions energy distribuion profile at the electrode obtained from 
PIC/MC simulator after running the simulator for 3000 periods and averaging 
over 150 periods. 

Fig. 7 (d) - Ions energy distribuion profile at the electrode obtained from 
PIC/MC simulator after running the simulator for 3000 periods and averaging 
over 600 periods. 
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Number of MC particles Number of MC particles 

50 100 150 50 100 150 200 250 



Ion density (1=1.0E8/CC) Ion density (1=1.0E8/cc) Ion density (1=1.0E8/cc) 

0.5 1.0 1.5 0.001 0.010 0.100 1.000 10.000 0.001 0.010 0.100 1.000 
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Electric field (V/cm) 
-100 -50 50 100 




Electric field (V/cm) 
-50 50 




Negative ion current density (A/m A 2) Ion current density (A/m A 2) 

-0.015 -0.005 0.005 0.015 -0.3 -0.1 0.1 0.2 0.3 
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